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Abstract 

In this paper, we consider the nonlinear dynamical behaviors of some tabu leaning neuron models. We 
first consider a tabu learning single neuron model. By choosing the memory decay rate as a bifurcation 
parameter, we prove that Hopf bifurcation occurs in the neuron. The stability of the bifurcating periodic 
solutions and the direction of the Hopf bifurcation are determined by applying the normal form theory. 
We give a numerical example to verify the theoretical analysis. Then, we demonstrate the chaotic 
behavior in such a neuron with sinusoidal external input, via computer simulations. Finally, we study 
the chaotic behaviors in tabu learning two-neuron models, with linear and quadratic proximity functions 
respectively. 
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1 Introduction 

Neurons as the fundamental elements of the brain can generate complex dynamical behaviors. The studies 
of nonlinear dynamics in both biological and artificial neural network models are very important in that, on 
the one hand, the studies may provide explanations for the rich dynamics including oscillations and chaos 
in biological neural networks, and on the other hand, if we understand more about the dynamical behaviors 
of various neural networks then we can design some control methods to achieve desirable system behaviors 
in artificial neural networks. 

In recent years, dynamical characteristics of various neural networks have become a focal subject of 
intensive research studies. Bifurcation and chaos have been identified and investigated in many kinds of 
neural networks. For example, chaotic solutions were obtained in [Kurten & Clark, 1986], from a neural 
network consisting of 26 neurons. In [Babcock & Westervelt, 1987], a two-neuron model was studied, with or 
without time delays. Numerical solutions of differential equations with electronic circuit models of chaotic 
neural networks were qualitatively compared in [Kepler et ah, 1990]. In [Das II et al., 1991], a chaotic 
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neural network with four neurons was investigated, and chaotic behavior was found in [Zou & Nossek, 1993] 
from a cellular neural network with three cells. In [Li et al., 2001], we also presented some findings of 
chaotic phenomenon in a three-neuron hysteretic Hopfield-type neural network. In [Bondarenko, 1997], a 
high-dimensional chaotic neural network under external sinusoidal force was studied. In [Ueta & Chen, 
2001], bifurcation and chaos as well as their control in a system of strongly connected neural oscillators were 
discussed. In [Das et al., 2002], chaos in a three-dimensional general neural network model was investigated. 
In [Olien & Bclair, 1997; Wei & Ruan, 1999; Gopalsamy et al., 1998; Liao et al., 1999; Liao et al., 2001; Liao 
et al., 2001; Gilli, 1993; Lu, 2002], the bifurcation and chaotic behaviors of various delayed neural networks 
were studied. Moreover, the chaotic phenomenon in a neural-network-based learning algorithm was reported 
in [Li et al., 2003]. 

Tabu learning [Beyer & Ogicr, 1991] applies the concept of tabu search [Glover, 1989; Glover, 1990] 
to neural networks for solving optimization problems. By continuously increasing the energy surface in a 
neighborhood of the current network state, it penalizes those states that have already been visited. This 
enables the state trajectory to climb out of local minima while tending toward those areas that have not 
yet been visited, thus performing an efficient search through the problem's solution space. Note that, unlike 
most existing neural-network-based methods for optimization, the goal of the tabu learning method is not 
to force the network state to converge to an optimal or nearly optimal solution, but rather, the network 
conducts an efficient search through the solution space. Knowing this, it is very natural for one to ask what 
the state trajectory of the tabu learning neural network is like? In the present paper, the aim is to provide 
an answer to this question by studying the dynamical behaviors of the tabu learning neural network. 

More precisely, this paper investigates the nonlinear dynamical behaviors of a tabu leaning single neuron 
model and of two two- neuron models, respectively. By choosing the memory decay rate as a bifurcation 
parameter, it is proved that Hopf bifurcation occurs in the tabu leaning single neuron model. The sta- 
bility of the bifurcating periodic solutions and the direction of the Hopf bifurcation are then determined 
by applying the normal form theory. Chaotic behavior of such a neuron with small sinusoidal external 
input is also demonstrated by computer simulations. To that end, the chaotic behaviors in tabu learning 
two- neuron models with linear and quadratic proximity functions, respectively, are studied by means of 
computer simulations. 

The paper is organized as follows. The tabu learning neural network is first described in Section 2. In 
Section 3, the Hopf bifurcation of a single neuron model is studied. In Section 4, the observed chaotic 
behavior of the single neuron model is described. In Section 5, the dynamical behaviors of a two-neuron 
model with a linear proximity function are studied. And, in Section 6, the dynamical behaviors of a two- 
neuron model with a quadratic proximity function are investigated. Finally, conclusions are drawn in Section 
7. 
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2 Tabu learning neural network 

The common problem of minimizing an objective function in many optimization schemes by using the 
Hopfield-type neural networks of the form 

C l ii l = —m + V TijVj + h (1) 

J 

can be translated into minimizing the following energy function: 

E ° = -\ E E WiVj - E 7 * y * + E ^ P / _1 ( s ) ds ( 2 ) 

i j i i 1 " 

where Vi = /(•) is the activation function, is the state of neuron i, d and i?i are positive constants, 

Tij represents the strength of the connection from neuron j to neuron i, and 7j represents the input current 
to neuron i. 

In the tabu leaning method, the energy surface E n is continuously increasing in a neighborhood of the 
current state. At time t, the energy surface is given by 

E t =E + F t (V) (3) 

with 



F t (V)=p[ e a(s - t) P(V,F(s))ds 
Jo 



where a and /? are positive constants and P(V, W) is a measure of the proximity of the two vectors V and 
W (thus, P(V, W) is maximized when V = W). Therefore, states "nearest" those already being visited are 
penalized the most, so that the penalty encourages a search through states that have not yet been visited. 
The exponential kernel keeps the integral from increasing to infinity and results in a higher penalty for 
states being visited most recently. The latter property also helps the network in quickly climbing out of local 
minima. 

In solving optimization problems, the memory decay rate a and the learning rate (3 must be chosen 
carefully. If a is too large, states already being visited are likely to be revisited again; but if it is too small, 
the network may require a long time to climb out of local minima. In addition, (3 must be chosen so as 
to achieve a balance between searching for the states that minimize Eq and searching for the states that 
minimize F t (V). 

One possible choice for the proximity function P(V, W) is defined by [Beyer & Ogier, 1991] 

Pi(V,W0 = ^(l + Wi) (4) 

i 

where Vi and Wi are the components of vectors V and W, respectively. This is called the linear proximity 
function for it is linear in V. 

A weakness of Pi(V, W) is that it penalizes too heavily those vectors that are not close to W. A proximity 
function, which is better in this respect, is given by [Beyer & Ogier, 1991] 

w w) = \ E(* + wx 1 + ( 5 ) 
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This is called the quadratic proximity function for it is quadratic in V, which leads to a quadratic penalty 
term F t (V). 

If the linear proximity function is selected, then the neural network that performs gradient descent on E t 
has the following state equation [Beyer & Ogier, 1991]: 



Cii- — -v • + V T V I I dFtp-0 

= -jrUi + YtiTijVj+Ii + Jiit) 



(6) 



where 



Mt) := -13 [ e a ^Vi(s)ds 
Jo 

Therefore, Ji satisfy the following learning equation: 

Ji = -aJi - (3Vi 

In the next two sections, a tabu learning single neuron model is studied, which is described by 

Cx = -x/R + af(x) +y + I 
y = -ay - f3f(x) 

where C, R, a, a and (5 are all positive constants, and / is the external input. 



(7) 



(8) 



3 Hopf bifurcation in a tabu learning single neuron model 

3.1 Existence of Hopf bifurcation 

In this subsection, conditions for the existence of a Hopf bifurcation in a tabu learning single neuron model 
is derived. The model is obtained by setting C = 1, R = 1, and / = in system (8), as follows: 



x = —x + af(x) + y 

y = -ay - Pf(x) 



(9) 



Suppose that in (9), / e C 3 (R), /(0) = 0, and that a, a and (3 satisfy the inequalities a > 0, b > 







X 




a/'(0) - 1 1 




X 


+ 


±af"(0)x 2 + ±af"{0)x 3 + --- 


y 




-(3f'(Q) -a 




y 







M < 1, where |/'(0)| < M. It is clear that system (9) have an unique equilibrium (0,0) under 
these conditions. 

Expanding system (9) into first, second, third and other higher-order terms about the equilibrium (0,0), 
one has 

(10) 



The associated characteristic equation of its linearized system is 

A 2 + (1 - o/'(0) + a)\ + (1 - o/'(0))a + f3f'(0) = 

Define 

61 =61 (a) = (l-a/'(0) + a) 
b 2 = 62(a) = (1 - o/'(0))a + (3f'(0) > 
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Then, the Routh-Hurwitz criterion implies that the equilibrium (0,0) of system (9) is locally asymptotically 
stable if 61(a) > 0. If 

a = af'(0) - 1 

then 61 («o) = 0, and the characteristic equation has one pair of purely imaginary roots, A^ 2 = icooh where 
= V^2(ao) = \J /3/'(0) — Oq. It follows from simple calculation that 



d[Re{Xi)] 



da 



«0 



1 

— 2<° 



The above analysis is summarized as follows: 

THEOREM 1: If a > a/'(0) - 1, then the equilibrium (0,0) of system (9) is locally asymptotically 
stable. If /3/'(0) — a 2 , > 0, then as a passes through the critical value a = af'(0) — 1, there is a Hopf 
bifurcation at the equilibrium (0,0). 



3.2 Stability of bifurcating periodic solutions 

In this subsection, the stability of the bifurcating periodic solutions is studied. The method used is based 
on the normal form theory [Hassard et al., 1981]. For notational convenience, let the above a = «o + 7, so 
7 = is the Hopf bifurcation value for system (9). 
The eigenvactor v\ associated with A x = + iuo is 



Vl 



H + iu- af (0) + 1 



Define 



and 



P = (Re«i, -Im«i) = 



1 

M-a/'(0) + l -u 



2/1 




X 


. V2 . 




. y . 



Then, the system of y\, yi is obtained from (10) as 

Vl = m - lo V2 + ia/"(0)y2 + ^af"(0)yf + ■■■ 

y 2 = um + mn + o(M " o/ ^ )+1)+/? r(0)i/? + a ^- af '£ )+1)+p r^yf + ■■■ 

Since one only needs to evaluate \ii , Ti , and P2 , which are obtained from Ci(0) alone [Hassard et al., 1981], 
one may set 

a = a Q = af'(0) - 1 

so that, in the above system of y\ and j/2, 



Thus, the system becomes 



0, and lu = lu = VWI - (a/'(0) - l) 2 

Vi = -u) V2 + F 1 (y 1 ,y 2 ;0) 
ill =u yi +i ? 2(2/i,J/2;0) 
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where 



and 



y 2 ;0) = \af"{Q)yl + \af"\<i)fi + 



F 2 {yi,y 2 ;0) = 



a ( i - afm +I + an- arm +p rm + . . 



2^o 



Next, the following quantities are evaluated at 
5ii = 



6uj 

(i.e. a = a ) and (2/1,1/2) = (0, 0). One has 



1 §111 4. j 

9yf ^ 9S/| ^ 



"s 2 f 2 , 


9 2 F 2 ~ 




. 9y? i 


9yl _ 
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502 = 



/"(0) 



a + z 



a(l-a/ , (0))+ff 



wo 

2 Fi _ o d 2 F 



|/"(0) 



520 



d 2 Fi 



dyidy 2 
a(l-af'(0))+l3 

a 2 F t _|_ 2 a 2 F 2 



+ i 



9 2 F 2 _ *9 2 F 2 



2 Fi 



a + z 

!)y. 



= l/"(0) 



a + z 



521 



dyidy 2 
a(l-a/ , (0))+/j 

3 3 Fi , 3 F 2 



+ i 



2 F 2 



dy| dyidy 2 



d 2 F 2 _ o Q 2 Fi 



(11) 



dyiSyJ dy'{dy 2 



d 6 F 2 

dy'i 



+ i 



5 F 2 



d 3 F 2 

dyidy'i 



d 3 F! 
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d 3 Fi 
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8 



f"(0) 



a + z 



,.a(l-n/'(0)) +) 3 



Based on the above analysis, one can see that each gij in (11) is determined by the parameters in (10). 
Thus, one can compute the following quantities: 



<)2L 



(12) 



Ci(0)= 2 ^[5205n-2| 5ll | 2 -i|502| 2 ] + 
fi 2 = -RcCi(0)//i'(0) = 2RcCi(0) 
r 2 = -(ImCi(0)+M2w'(0))/wo 
#2 = 2RcCi(0) 

Now, the main results of this subsection are summarized as follows. 

THEOREM 2: In formulas (12), fi 2 determines the direction of the Hopf bifurcation: if fi 2 > (< 0), 
then the Hopf bifurcation is supercritical (subcritical) and the bifurcating periodic solutions exist for a > 
ctQ (< ao); P2 determines the stability of the bifurcating periodic solutions: the solutions are orbitally stable 
(unstable) if /3 2 < (> 0); and t 2 determines the period of the bifurcating periodic solutions: the period 
increases (decreases) if r 2 > (< 0). 



3.3 A numerical example 

In this subsection, an example in the form of system (9) is discussed, with a = 1.6, /(•) = tanh(-) and (3 = 0.5. 
By Theorem 1, one can determine that 

a Q = 0.6, loq = 0.3742 

It follows from the results in Subsection 3.2 that 

H2 = -0.4, r 2 = 0.9446, (3 2 = -0.4 

These calculations prove that the equilibrium (0, 0) is stable when a > a>o, as shown by Fig. 1, where a = 0.7. 
When a passes through the critical value a = 0.6, the equilibrium losses its stability and a Hopf bifurcation 
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occurs, i.e., a family of periodic solutions bifurcate out of the equilibrium. Each individual periodic orbit is 
stable since < 0. Since fi 2 < 0, the bifurcating periodic solutions exist at least for values of a slightly less 
than the critical value. Choosing a = 0.02, as predicted by the theory, Fig. 2 shows that there indeed is a 
stable limit cycle. Since r 2 > 0, the period of the periodic solutions increases as a increases. For a = 0.5, 
the phase plot and the waveform plot are shown in Fig. 3. Comparing Fig. 2 with Fig. 3, one can conclude 
that the period of a = 0.5 is longer than that of a = 0.02. 




100 



Figure 1: Phase plot and waveform plot for system (9) with a = 0.7 




Figure 2: Phase plot and waveform plot for system (9) with a = 0.02 



4 Chaos in an excited tabu learning single neuron model 

In this section, the dynamical behaviors of a tabu learning single neuron model, with external sinusoidal 
input, are studied. 

In system (8), let I be 

I(t) = esin(ujt) 

and let the nonlinear activation function be 

f{x)=pxe- (px ^^ 2 (13) 

where e, u), p, and a are positive constants. Also, in system (8), set the parameters be C = 1, i? = 7, a = 
0.3, e = 0.2, uo = 27r, a = 0.6, (i = 0.9, and set p = 8, a 2 = 0.2 in (13). It should be noted that several other 
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Figure 3: Phase plot and waveform plot for system (9) with a = 0.5 
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Figure 4: Phase plot 



Figure 5: Waveform plot of x 



parameters have also been examined, showing similar dynamical phenomena. Due to the space limitation, 
those results are not presented here. 

The phase plot of x and y and the waveform plot of x are shown in Fig. 4 and Fig. 5, respectively. And 
in Fig. 6, the power spectrum of the neuron state x, calculated using the general FFT, is plotted. 

In order to calculate the largest Lyapunov exponent, the method introduced in [Wolf et al., 1985] is used, 
in which one monitors the long-term evolution of a single pair of nearby orbits. The largest Lyapunov 
exponent is defined as 



A 



1 M 



tM — to 



fe=l 



L'(t k ) 
L(t 



k—l ) 



where M is the total number of replacement steps, L(tk—i) is the distance between the two initial points at 
the time instant tk-i- After a time step S — tk — ifc-i, the initial length will have evolved to length L'(tk), 
as detailed in [Wolf et al., 1985]. In the simulations, the largest Lyapunov exponent A is calculated from a 
time series of A = 100, 000 points. The largest Lyapunov exponent is obtained as A = 0.7030. 

From Figs. 4-6 and the largest Lyapunov exponent, one can see that there exists chaotic behavior in this 
nonautonomous tabu learning single neuron model. 
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Figure 6: Power spectrum of the time series of x 

5 Chaos in a two-neuron model with a linear proximity function 

Now, consider a two-neuron model with a linear proximity function. From (6) and (7), one knows that there 
are totally four differential equations in this system: 



C\ii\= -^u 1 +T 11 V 1 +T 12 V 2 + h+ Ji 

C 2 U2= -^u 2 + T 21 V 1 +T 22 V 2 + I 2 + h 

Ji = -aJi - 0Vi 

J 2 = -aJ 2 - (3V 2 



(14) 



In the following simulations, without loss of generality, let C\ — C 2 — 1, R\ = R 2 = 10, h = I 2 = 0, 
a = 0.1, (3 be a variable parameter, and let the weight matrix be 



T = 



0.1 0.5 

-1 2 



Moreover, let the activation function be Vi — f(ut) = tanh(5ui). It should be noted that several other 
parameters have also been examined, showing similar dynamical behaviors. Due to space limitation, those 
results are not presented here. 
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Figure 7: Phase plot and waveform plot of the two-neuron model (14) (/? = 0.1) 
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Figure 8: Phase plot and waveform plot of the two-neuron model (14) {(3 = 0.5) 
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Figure 9: 3-D and 2-D phase plots of the two-neuron model (14) {(3 = 1) 

Start from the case of (3 = 0.1, with which the neural network produces periodic solutions. Fig. 7 shows 
the 3-D phase plot and the waveform plot of u\. Fig. 8 shows the 3-D phase plot and the waveform plot of 
the network orbit when (3 — 0.5. In this case, it is also a periodic solution, but it is different from that of the 
case (3 = 0.1. When (3=1, the orbit of the neural network is chaotic, and the largest Lyapunov exponent is 
0.2419. The 3-D and 2-D plots of the chaotic attractor are shown in Fig. 9. The corresponding waveform 
and power spectrum are shown in Fig. 10, where the spectrum has been truncated at ^=0.5 HZ. 



6 Chaos in a two-neuron model with a quadratic proximity func- 
tion 

If the quadratic proximity function (5) is selected, then the tabu learning neural network that performs 
gradient descent minimization on E t has the following state equations: 

Cu - 1 u i y T V I I dFt(y) 

= + + SijityVj +I t + Ji{t) 
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Figure 10: Waveform plot and power spectrum of the two-neuron model (14) (fi = 1) 



where 



Su(t) 
Ji(t) 



:= 

:= -0(n-l)S*e a ^Vi(8)d8 
in which n is the number of neurons. Therefore, Sij(t), i ^ j, and Ji(t) satisfy the following learning 
equations: 

Sij = -aSij - pViVj (16) 



(18) 



Ji = -aJ l - (3{n - l)Vi (17) 
There are totally six differential equations in this model: 

Ciui = -±-u x + (Tia + S 12 (t))V 2 +h + ,h{t) 
C 2 u 2 = -^u 2 + (T 2 i + S 2 i(t))V 1 +I 2 + J 2 (t) 
S12 = -0S12 - PViV 2 
S21 = -0S21 - /3V 2 Vi 
Ji = -aJi - /3Fi 
j 2 = -a J2 - /8V2 

For this model, only its chaotic behavior is discussed. Let C\ = C 2 = 1, R\ = R 2 = 10, Ji = J2 = 0, a = 0.1, 
/3 = 100, Vi = tanh(lOuj), and the weight matrix be 

30 
50 

Then, the neural network is chaotic. The 3-D and 2-D phase plots of the chaotic attractor are shown in Fig. 
11 and Fig. 12, respectively. The waveform plot and power spectrum of u\ are shown in Fig. 13, in which 
the largest Lyapunov exponent is 0.3160. 
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Figure 11: 3-D phase plot of the chaotic attractor of the two-neuron model (18) 

7 Conclusions 

The nonlinear dynamical behaviors of some tabu learning neuron models have been studied in this paper. 
By choosing the memory decay rate as a bifurcation parameter, it has been proved that Hopf bifurcation 
occurs in the single neuron model. The stability of bifurcating periodic solutions and the direction of the 
Hopf bifurcation are determined based on the normal form theory. From the waveform diagrams, the phase 
plots, the power spectra, and the largest Lyapunov exponents, one can find chaotic phenomena in the single 
neuron model with small external sinusoidal input and in the two-neuron models. These tabu learning 
models, although simple, have rich complex dynamics, and deserve further investigations. 
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